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Abstract. We present Lie group integrators for nonlinear stochastic differential equations with 
non-commutative vector fields whose solution evolves on a smooth finite dimensional manifold. Given 
a Lie group action that generates transport along the manifold, we pull back the stochastic flow on 
the manifold to the Lie group via the action, and subsequently pull back the flow to the corresponding 
Lie algebra via the exponential map. We construct an approximation to the stochastic flow in the 
Lie algebra via closed operations and then push back to the Lie group and then to the manifold, 
thus ensuring our approximation lies in the manifold. We call such schemes stochastic Munthe-Kaas 
methods after their deterministic counterparts. We also present stochastic Lie group integration 
schemes based on Castell— Gaines methods. These involve using an underlying ordinary differential 
integrator to approximate the flow generated by a truncated stochastic exponential Lie series. They 
become stochastic Lie group integrator schemes if we use Munthe-Kaas methods as the underlying 
ordinary differential integrator. Further, we show that some Castell-Gaines methods are uniformly 
more accurate than the corresponding stochastic Taylor schemes. Lastly we demonstrate our methods 
by simulating the dynamics of a free rigid body such as a satellite and an autonomous underwater 
vehicle both perturbed by two independent multiplicative stochastic noise processes. 
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1. Introduction. We are interested in designing Lie group numerical schemes 
for the strong approximation of nonlinear Stratonovich stochastic differential equa- 
tions of the form 



Here W 1 , . . . , W are d independent scalar Wiener processes and W° = t. We suppose 
that the solution y evolves on a smooth n-dimensional submanifold M. of M. N with 
n < N and Vi : Ai x R + — > TAi, i — 0,1, ... ,d, are smooth vector fields which in 
local coordinates are Vi — J2j=i ^l^yy The flow-map <pt : M. — > M. of the integral 
equation is defined as the map taking the initial data yo to the solution y t at 

time t, i.e. y t = (fit ° yo- 

Our goal in this paper is to show how the Lie group integration methods developed 
by Munthe-Kaas and co-authors can be extended to stochastic differential equations 
on smooth manifolds (see Crouch and Grossman [5] and Munthe-Kaas [10]). Suppose 
we know that the exact solution of a given system of stochastic differential equations 
evolves on a smooth manifold M. (see Malliavin [35] or Emery [2]), but we can only 
find the solution pathwise numerically. How can we ensure that our approximate 
numerical solution also lies in the manifold? 

Suppose we are given a finite dimensional Lie group Q and Lie group action A yo 
that generates transport across the manifold M. from the starting point y a £ Ai via 
elements of Q . Then with any given elements £ in the Lie algebra g corresponding to 
the Lie group Q , we can associate the infinitesimal action using the Lie group action 
Aj, . The map £ ► is a Lie algebra homomorphism from g to X(A4), the Lie algebra 
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of vector fields over the manifold M.. Further the Lie subalgebra {A^ S X(A4) : £ € g} 
is isomorphic to a finite dimensional Lie algebra with the same structure constants 
(see Olver g2J, p. 56). 

Conversely, suppose we know that the Lie algebra generated by the set of govern- 
ing vector fields Vi, i = 0, 1, . . . , d, on M is finite dimensional, call this Xf(-M)- Then 
we know there exists a finite dimensional Lie group Q whose Lie algebra g has the 
same structure constants as Xf(M) relative to some basis, and there is a Lie group 
action A yo such that Vi = i = 0, 1, . . . , d, for some & £ g (see Olver [i2] . p. 56 or 
Kunita [3U], p. 194). The choice of group and action is not unique. 

In this paper we assume that there is a finite dimensional Lie group Q and action 
A Vo such that our set of governing vector fields Vi, i = 0,1, . . . ,d, are each infinitesimal 
Lie group actions generated by some element in g via A Bo , i.e. Vi = A^ for some £j 6 g, 
i = 0, 1, . . . , d. They are said to be fundamental vector fields. This means that we can 
write down the set of governing vector fields for a system of stochastic differential 
equations on the Lie group Q that, via the Lie group action A yo , generates the flow 
governed by the set of vector fields Vi on the manifold. The vector fields Vi on A4 are 
simply the push forward of the vector fields X^ on Q via the Lie group action Aj, . 
Typically the flow on the Lie group also needs to be computed numerically. We thus 
want the approximation to remain in the Lie group so that the Lie group action takes 
us back to the manifold. 

To achieve this, we pull back the set of governing vector fields X^ on Q to the set 
of governing vector fields on g, via the exponential map 'exp' from g to Q . Thus 
the stochastic flow generated on g by the vector fields generates the stochastic flow 
on Q generated by the X^ i . The set of governing vector fields on g are for each a G g: 

oo 
fc=0 

where is the fcth Bernoulli number and the adjoint operator ad CT is a closed operator 
on g, in fact ad CT o ( = [a, £], the Lie bracket on g. Now the essential point is that 
6 § and so the series on the right or any truncation of it is closed in g. Hence 
if we construct an approximation to our stochastic differential equation on g using 
the vector fields or an approximation of them achieved by truncating the series 
representation, then that approximation must reside in the Lie algebra g. We can then 
push the approximation in the Lie algebra forward onto the Lie group and then onto 
the manifold. Provided we compute the exponential map and action appropriately, our 
approximate solution lies in the manifold (to within machine accuracy) . In summary, 
for a given £ € g and any yo € M. we have the following commutative diagram: 

exp, v , n ^ ( A ™)* v /kA\ 

g ► X(G) ► X(M) 

«<T T*« j\ 

► y > M 

We have implicitly separated the governing set of vector fields Vi, i = 0, 1, . . . , d, 
from the driving path process w = (W 1 , . . . , W d ). Together they generate the unique 
solution process y £ M. to the stochastic differential equation ([LTJ) , When there 
is only one driving Wiener process (d = 1) the ltd map w <— > y is continuous in 
the topology of uniform convergence. When there are two or more driving processes 
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(d > 2) the Universal Limit Theorem tells us that the Ito map w <— > y is continuous 
in the p- variation topology, in particular for 2 < p < 3 (see Lyons [32] . Lyons and 
Qian [33] and Malliavin [35]). A Wiener path with d > 2 has finite p- variation for 
p > 2. This means that from a pathwise perspective, approximations to y constructed 
using successively refined approximations to w are only guaranteed to converge to 
the correct solution y, if we include information about the Levy chordal areas of 
the driving path process. Note however that the I/ 2 -norm of the 2-variation of a 
Wiener process is finite. In the Lie group integration procedure prescribed above we 
must solve a stochastic differential system on the Lie algebra g defined by the set 
of governing vector fields and the driving path process w = (W 1 , . . . , W ). In 
light of the Universal Limit Theorem and with stepsize adaptivity in mind in future 
(see Gaines and Lyons [20]), we for instance use in our examples order 1 stochastic 
numerical methods — that include the Levy chordal area — to solve for the flow on the 
Lie algebra g. 

We have thus explained the idea behind Munthc-Kaas methods and how they can 
be generalized to the stochastic setting. The first half of this paper formalizes this 
procedure. 

In the second half of this paper, we consider autonomous vector fields and con- 
struct stochastic Lie group integration schemes using Castell-Gaines methods. This 
approach proceeds as follows. We truncate the stochastic exponential Lie series expan- 
sion corresponding to the flow ip t of the solution process y to the stochastic differential 
equation We then approximate the driving path process w = (W 1 , . . . , W d ) by 

replacing it by a suitable nearby piecewise smooth path in the appropriate variation 
topology. An approximation to the solution y t requires the exponentiation of the 
approximate truncated exponential Lie series. This can be achieved by solving the 
system of ordinary differential equations driven by the vector field that is the approx- 
imate truncated exponential Lie series. If we use ordinary Munthe-Kaas methods as 
the underlying ordinary differential integrator the Castell-Gaines method becomes a 
stochastic Lie group integrator. 

Further, based on the Castell-Gaines approach we then present uniformly accurate 
exponential Lie series integrators that are globally more accurate than their stochastic 
Taylor counterpart schemes (these are investigated in detail in Lord, Malham and 
Wiese [3JJ for linear stochastic differential equations). They require the assumption 
that a sufficiently accurate underlying ordinary differential integrator is used; that 
integrator could for example be an ordinary Lie group Munthe-Kaas method. In 
the case of two driving Wiener processes we derive the order 1/2, and in the case 
of one driving Wiener process the order 1 uniformly accurate exponential Lie series 
integrators. As a consequence we confirm the asymptotic efficiency properties for 
both these schemes proved by Castell and Gaines [5] (see Newton [UJ for more details 
on the concept of asymptotic efficiency). We also present in the case of one driving 
Wiener process a new order 3/2 uniformly accurate exponential Lie series integrator 
(also see Lord, Malham and Wiese [3T]). 

We present two physical applications that demonstrate the advantage of using 
stochastic Munthe-Kaas methods. First we consider a free rigid body which for ex- 
ample could model the dynamics of a satellite. We suppose that it is perturbed by 
two independent multiplicative stochastic noise processes. The governing vector fields 
are non-commutative and the corresponding exact stochastic flow evolves on the unit 
sphere. We show that the stochastic Munthe-Kaas method, with an order 1 stochastic 
Taylor integrator used to progress along the corresponding Lie algebra, preserves the 



4 



Malham and Wiese 



approximate solution in the unit sphere manifold to within machine error. However 
when an order 1 stochastic Taylor integrator is used directly, the solution leaves the 
unit sphere. The contrast between these two methods is more emphatically demon- 
strated in our second application. Here we consider an autonomous underwater vehicle 
that is also perturbed by two independent multiplicative stochastic noise processes. 
The exact stochastic flow evolves on the manifold which is the dual of the Euclidean 
Lie algebra se(3); two independent Casimirs are conserved by the exact flow. Again 
the stochastic Munthe-Kass method preserves the Casimirs to within machine error. 
However the order 1 stochastic Taylor integrator is not only unstable for large step- 
sizes, but the approximation drifts off the manifold and makes a dramatic excursion 
off to infinity in the embedding space ]R 6 . 

Preserving the approximate flow on the manifold of the exact dynamics may be a 
required property for physical or financial systems driven by smooth or rough paths — 
for general references see Iserles, Munthe-Kaas, N0rsett and Zanna [25], Hairer, Lubich 
and Wanner [22] , Elworthy [13] , Lyons and Qian [33] and Milstein and Tretyakov [38 . 
Stochastic Lie group integrators in the form of Magnus integrators for linear stochastic 
differential equations were investigated by Burrage and Burrage [5]. They were also 
used in the guise of Mobius schemes (see Schiff and Shnider 03]) to solve stochastic 
Riccati equations by Lord, Malham and Wiese [31] where they outperformed direct 
stochastic Taylor methods. Further applications where they might be applied include: 
backward stochastic Riccati equations arising in optimal stochastic linear-quadratic 
control (Kohlmann and Tang 28 1); jump diffusion processes on matrix Lie groups 
for Bayesian inference (Srivastava, Miller and Grenander [H]); fractional Brownian 
motions on Lie groups (Baudoin and Coutin [3]) and stochastic dynamics triggered 
by DNA damage (Chickarmane, Ray, Sauro and Nadim [TO]). 

Our paper is outlined as follows. In Section 2 we present the basic geometric setup, 
sans stochasticity. In particular we present a generalized right translation vector field 
on a Lie group that forms the basis of our subsequent transformation from the Lie 
group to the manifold. Using a Lie group action, this vector field pushes forward 
to an infinitesimal Lie group action vector field that generates a flow on the smooth 
manifold. In Section 3 we specialize to the case of a matrix Lie group and using 
the exponential map, derive the pullback of the generalized right translation vector 
field on the Lie group to the corresponding vector field on the Lie algebra. To help 
give some context to our overall scheme, we provide in Section 4 illustrative examples 
of manifolds and natural choices for associated Lie groups and actions that generate 
flows on those manifolds. Then in Section 5 we show how a flow on a smooth manifold 
corresponding to a stochastic differential equation can be generated by a stochastic 
flow on a Lie algebra via a Lie algebra action. We explicitly present stochastic Munthe- 
Kaas Lie group integration methods in Section 6. We start the second half of our 
paper by reviewing the exponential Lie series for stochastic differential equations in 
Section 7. We show in Section 8 how to construct geometric stochastic Castell-Gaines 
numerical methods. In particular we also present uniformly accurate exponential Lie 
series numerical schemes that not only can be used as geometric stochastic integrators, 
but also are always more accurate than stochastic Taylor numerical schemes of the 
corresponding order. In Section 9 we present our concrete numerical examples. Finally 
in Section 10 we conclude and present some further future applications and directions. 

2. Lie group actions. Suppose M. is a smooth finite n-dimensional submanifold 
of 1 N with n < N. We use X(A4) to denote the Lie algebra of vector fields on the 
manifold A4, equipped with the Lie-Jacobi bracket [{/, V] = U ■ W — V ■ VJ7, for all 
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U, V £ X(Ai). Let ^ denote a finite dimensional Lie group. 

Definition 2.1 (Lie group action). A left Lie group action of a Lie group Q on a 
manifold Ai is a smooth map A: Q x Ai — > satisfying for all y £ Ai and R, S £ Q: 
(1) A(id,y) = y; (2j A{R,A{S,y)) = A{RS,y). We denote A y oS = A(S,y). 

Hereafter we suppose yo £ is fixed and focus on the action map A yo : <?— >A4. 
We assume that the Lie group action A is transitive, i.e. transport across the manifold 
from any point yo £ Ai to any other point y £ Ai can always be achieved via a group 
element S £ Q with y = A yo o S (Marsden and Ratiu [37], p. 310). 

We define the Lie algebra g associated with the Lie group Q to be the vector space 
of all right invariant vector fields on Q . By standard construction this is isomorphic 
to the tangent space to Q at the identity id = \Ag (see Olver [32], p. 48 or Marsden 
and Ratiu [37, p. 269). 

Definition 2.2 (Generalized right translation vector field). Suppose we are 
given a smooth map £: Ai— >q. With each such map £ we associate a vector field 
AQ : Q — > X(G) defined as follows 

X^oS= d T exp(raA yo oS))S\ T=0 

for S £ Q, where 'exp' is the usual local diffeomorphism exp: q — > Q from a neigh- 
bourhood of the zero element o £ q to a neighbourhood of id £ Q . 

Definition 2.3 (Infinitesimal Lie group action). We associate with each vector 
field : G^X(G) a vector field : M— >X(M) as the push forward of X^ from Q to 
M. by A yo , i.e. = (Aj, )^A^, so that if S £ Q and y = A yo o S £ M., then 

h°V= d rA Va °7(t)I t=0 ' 

where j(t) £ Q , 7(0) = S and 9 T 7(r) = X^ o j(t) (the flow generated on Q by the 
vector field X^ starting at S £ Q). Naturally, as a vector field is linear, and also 

h°V = oA yo oS > 

the Lie derivative of A Vo along X^ at S £ Q. 
Remarks. 

1. The map A{S): Ai^Ai defined by y i— ► A(S) o y = A y o S represents a 
flow on A4. Hence if y = A(S) o y , the push forward of by A(S) is given by 
(A(S)),\£ = A A d s c (Marsden and Ratiu [37 , p. 317). 

2. We define the isotropy subgroup at yo £ Ai by Q yo = {S £ Q: A VQ oS = yo}; it 
is a closed subgroup of Q (see Helgason [53], p. 121 or Warner [35], p. 123). We define 
the global isotropy subgroup by Q M = r\ yo£M G yo = {S £ Q : A yo o S = y , Vj/ S M}\ 
it is a normal subgroup of Q (see Olver [35], p. 38). 

3. A Lie group action is said to be is effective /faithful if the map S <— > A(S) from 
Q to Diff(Ai), the group of diffeomorphisms on Ai, is one-to-one. This is equivalent 
to the condition that different group elements have different actions, i.e. Qm = {idg}. 
A Lie group action is said to be free if Q Vo = {idg} for all yo £ Ai, i.e. A Vo is a 
diffeomorphism from Q to Ai. For more details see Marsden and Ratiu [37], p. 310 
and Olver A2\, p. 38. 

4. The map 7: G/Gy —>A4 defined by 7: S ■ G yo *— * A yo o S is a diffeomorphism, 
i.e. Ai = G/ G Vo for any y £ Ai (a manifold Ai with a Lie group action A : G x .M-^Ai 
defined over it is thus diffeomorphic to a homogeneous manifold; see Warner [48], 
p. 123 or Olver [42], p. 40). Further, the induced action of G/Gm on -M is effective. 
Hence if A is not an effective action of G, we can replace it (without loss of generality) 
by the induced action of G/Gm ( see Olver [35], p. 38). 
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5. Our definition for the generalized right translation vector field X% on Q is 
motivated by the standard right translation vector field used to identify g, the vector 
space of right invariant vector fields on G, with T\^G, the tangent space to Q at the 
identity. When £ G g is constant, X^ G %(G) is right invariant and a Lie bracket 
on T\&Q can be defined via right extension by the corresponding Lie-Jacobi bracket 
for the vector fields X^ on X(Q). Unless £ G 9 is constant, X% is not in general 
right invariant. For further details see Varadarajan [47], Olver }42j . or Marsden and 
Ratiu [37]. 

6. The infinitesimal generator map ^ 1 — > from g to X(A4) is a Lie algebra 
homomorphism. If we identify g as the vector space of left invariant vector fields on Q 
this map becomes an anti-homomorphism. The Lie-Jacobi bracket as defined above 
gives the right (rather than left) Lie algebra stucture over the group of diffeomorphisms 
on M. . If in addition we take the Lie-Jacobi bracket to be minus that defined above — 
associated with the left Lie algebra structure — then the infinitesimal generator map 
becomes a homomorphism again. See for example Marsden and Ratiu 37J, p. 324 or 
Munthe-Kaas [40] . 

7. The image of g under the infinitesimal generator map £ 1— > forms a finite 
dimensional Lie algebra of vector fields on M. which is isomorphic to the Lie algebra of 
the effectively acting quotient group G/Gm ( see Olver [55], p. 56). Thus the tangent 
space to M. at any point is g and M inherents a connection from G /Gm ■ Connections 
are necessary to define martingales on manifolds, but not for defining semimartingales 
(our focus here); see Malliavin [36] and Emery [13] . 

8. A comprehensive study of the systematic construction of symmetry Lie groups 
from given vector fields can be found in Olver [42j . 

9. We assumed above that the vector fields Xf and \ are autonomous. However 
all results in this and subsequent sections up to Section [JJ can be straightforwardly 
extended to non-autonomous vector fields generated by £: M x R— >g with (y, t) 1— > 
f(y, t) for all y G M and tel. 

10. For full generality we want to suspend reference to embedding spaces as far as 
possible. However in subsequent sections to be concise we will more explicitly reclaim 
this context. 

3. Pull back to the Lie algebra. For ease of presentation, we will assume in 
this section that G is a matrix Lie group. Recall that the exponential map exp: g — ► G 
is a local diffeomorphism from a neighbourhood of o G g to a neighbourhood of id G G- 
Let ve : g— >g be the pull back of the vector field X^ : G~^%(G) from G to 9 via the 
exponential mapping exp: 9— >G, i.e. v^o a = exp* o a. If a G 9 then 



Here dexp~ : g—>g is the inverse of the right-trivialized tangent map of the exponential 
dexpg. : g^9 defined as follows. If /3(r) is a curve in 9 such that (3(0) — a and 
/3'(0) = r) G 9 then dexp: g x g^g is the local smooth map (Varadarajan [47], p. 108) 




(3.1) 



dexp^ o 77 = d T exp (3(t) \ T=0 
/ exp(ad CT ) — id 




exp(— a) 



Note that as a tangent map dexp CT : g-^g is linear. The inverse operator dexp^ is the 
operator series (jl.2p generated by considering the reciprocal of dexp CT . 
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To show that f)3 . 1 [) is true, if exp: g— *Q with iti->5 = exper, and (3(t) £ g with 
(3(0) = a and 9 T /3(r) = o /3(r), then: 

exp, v i oS= d T exp /3(r) | T=0 

= (dexp CT o v£ o a) exp(er) 

Since 'exp' is a diffcomorphism in a neighbourhood of o e g, this push forward calcu- 
lation establishes the pull back (|3.1[) for all a £ g in that neighbourhood. 

4. Illustrative examples. Suppose the vector field V : M. x R^X(.A4) gener- 
ates a flow solution y t £ M starting from y £ A4. Then assume there exists a: 

1. Lie group Q with corresponding Lie algebra g; 

2. Lie group action A ya : Q-^M. for which a starting point yo £ M. is fixed; 

3. Vector field : x R— >X(.A4) such that: V = A^, i.e. V is a fundamental 
vector field corresponding to the action A yo . 

Let us suppose Q is a matrix Lie group (or can be embedded into a matrix Lie 
group, for example the Euclidean group SE(3) is naturally embedded into the special 
linear group SX(4;R)). We have for all S £ G and t € R, 

X ( (S,t)^t(A yo (S),t)S. (4.1) 

If V = Ac for some £ : M— >g, some Lie group Q and corresponding action A yo , then 
the flow generated by X^ on Q drives the flow generated by F on A4. In each of the 
examples below, given the manifold A4, we present a natural Lie group and action 
associated with the manifold structure, and identify vector fields which generate flows 
on the manifold via the Lie group. 

Stiefel manifold V n , fe . Suppose M = V„, fe = {y £ R nxk : y T y = I}. Take 
Q = SO(n), the special orthogonal group, and A yo (S) = Syo, the action of left 
multiplication. The corresponding Lie algebra g = so(n). Then by direct calculation 
\(y) — £(y>t)y- Hence if the given vector field V(y,t) = £(y,t)y, then the push 
forward of the flow generated by X$ (S, t) on Q in (|4.ip is the flow generated by V on 
M. Note that the unit sphere § 2 = ¥3,1, i.e. § 2 is just a particular Stiefel manifold. 
In Section 9 as an application, we consider rigid body dynamics evolving on S 2 . 

Isospectral manifold §„. Suppose M. = § n = {y G R nxn : y T = y}, the set of 
n x n real symmetric matrices. Take Q = 0(n), the orthogonal group and A yo (S) = 
SyaS T , which is an isospectral action (Munthe-Kaas [10]). The corresponding Lie 
algebra is g = so(n). Again, by direct calculation X^(y) = £(y,t)y — y£(y,t). Hence 
if the given vector field V (y, t) = £ (y, t)y — y(;(y,t), then the push forward of the flow 
generated by X$(S, t) on Q in (14. ip is the flow generated by V on M.. 

Dual of the Euclidean algebra se(3)*. Suppose M — se(3)* = R 3 , the dual 
of the Euclidean algebra se(3) of the Euclidean group SE(3) = {(s,p) £ SE(3) : s £ 
50(3), p S R 3 }. Take Q = SE(3) so g = se(3) and A = Ad*: Q x g*^g*, the 
coadjoint action of Q on g*. Then by direct calculation A^(y) — — ad^(y). Since 
A^(y) in linear in £ and — Ae(y) = A_^(y), it follows that if V(y) — ad*-(y), then 
the push forward of the flow generated by X-^(S, t) = —£(A yo (S),t) S on Q is the 
flow generated by V on M.. For more details see Section 9 where we investigate the 
dynamics of an autonomous underwater vehicle evolving on se(3)*. 
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Grassmannian manifold Gr(fc, n). The Grassmannian manifold M. = Gr(fc, n) 
is the space of /c-dimensional subspaces of M. n . Take Q = GL(n), the general linear 
matrix group, where if S € GL(n), we identify 

a (3 
7 S 

where the block matrices a, (3, 7 and 5 are sizes k x k, k x (n — k), (n — k) x k and 
(n — k) x (n — fc), respectively (see Schiff and Shnider [43]; Munthe-Kaas [40]). We 
choose the action of GL(n) on Gr(fc,n) to be the generalized Mobius transformation 
A yo (S) = (ay Q + 0)(<yy o + 5)- 1 . Hence if 

b (t)\ 

W [c(t) d(t)J ' 

then direct calculation reveals that A^(y) = a(t)y + b(t) — yc{t)y — yd{t). Hence if the 
given vector field V(y) — a(t)y + b(t) — yc(t)y — yd(t), then the push forward of the 
flow generated by X^(S, t) — £(t) S on Q is the flow generated by V on Gr(fc, n). 

5. Stochastic Lie group integration. We show that if a Lie group action 
A : Q x A4— >M. exists, then for yo € M. fixed, the Lie algebra action A yo o exp : g— >M 
carries a flow on g to a flow on M.. 

Theorem 5.1. Suppose there exists a Lie group action A: Q x .M— >A4. Then if 
there exists a process a € g and a stopping time T* such that on [0, T*), a satisfies 
the Stratonovich stochastic differential equation 

d r-t 

o-t=J2 ifcooydW?, (5-1) 

then the process y = A yo o exp a ^ M. satisfies the Stratonovich stochastic differential 
equation on [0, T*): 

Vt = Vo + Y, f\<oy T dWi. (5-2) 
Proof. Using Ito's lemma, if <r t S g satisfies (|5.1[) then A ao o exp<r t satisfies 

d r-t 

A yo o exp CT t = A yo o exp o + V" / £„ e . o A yo o exp cr T dW* . 



via 



Now recall that for each i = 0, 1, . . . , d, is the push forward of from g to Q . 
the exponential map, and that A^ is the push forward of X^ from Q to .M via A 
and so the Lie derivative 

A) €< A ao o exp o- t ee A 5i o y t . 

Then since yt = A yo o exp at , we conclude that y e .M is a process satisfying the 
stochastic differential equation (|5.2p . 

Corollary 5.2. Suppose that for each i = 0, 1, . . . , d there exists £j : M— >q such 
that the vector field Vi : M.— >X(M) and A^ : SA—*X(JA) can be identified, i.e. 

Vi ee A c , . (5.3) 
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Then the push forward by 'h yo oexp ' of the flow on the Lie algebra manifold q generated 
by the stochastic differential equation (15. ip is the flow on the smooth manifold M. 
generated by the stochastic differential equation (|5.2p . whose solution can be expressed 
in the form yt = A yo o exp a t . 

Remark. If the action is free then 'A yo o exp' is a diffeomorphism from a neigh- 
bourhood of o £ 9 to a neighbourhood of yo £ A4. 

6. Stochastic Munthe-Kaas methods. Assuming that the vector fields in our 
original stochastic differential equation (jl.ip are fundamental and satisfy (|5.3[) . then 
stochastic Munthe-Kaas methods are constructed as follows: 

1. Subdivide the global interval of integration [0,T] into subintervals [i n) i n +i]. 

2. Starting with to = 0, repeat the next two steps over successive intervals 
[tm tn+i] until t n+ i = T. 

3. Compute an approximate solution &t n j n+1 to (|5.ip across [t n ,t n+ i] using a 
stochastic Taylor, stochastic Runge-Kutta or Castell-Gaines method. 

4. Compute the approximate solution yt n+1 ~ ^ Vtn ° expa tn ,t n+1 - 

Note that by construction at n ,t n+1 & 9 because the stochastic differential equa- 
tion (|5.ip (or any stochastic Taylor or other sensible approximation) evolves the so- 
lution locally on the Lie algebra g via the vector fields : 0—>0. Suitable methods 
for approximating the exponential map to ensure it maps g to Q appropriately can be 
found in Iserles and Zanna |26J. Then by construction yt n+1 € M- 

For example, with two Wiener processes and autonomous vector fields o a, an 
order 1 stochastic Taylor Munthe-Kaas method is based on 

<?t„,t„+i = (Jqv^ + Jiv^+ J2V^ + \Jlv\ x + Ji2V£ 1 V£ 2 + J2iv^ 2 v^ + \jlvl 2 )oo , (6.1) 

evaluated at the zero element 065. Typically 'dexp^ 1 ' is truncated to only include 
the necessary low order terms to maintain the order of the numerical scheme. 

Remark. It is natural to invoke Ado's Theorem (see for example Olver [42] p. 54): 
any finite dimensional Lie algebra is isomorphic to a Lie subalgebra of flt(n) (the 
general linear algebra) for some n £ N. However as Munthe-Kaas [40] points out, 
directly using a matrix representation for the given Lie group might not lead to the 
optimal computational implementation (other data structures might do so). 

7. Exponential Lie series. The stochastic Taylor series is known in different 
contexts as the Neumann series, Peano-Baker series or Feynman-Dyson path ordered 
exponential. If the vector fields in the stochastic differential equation (jl.ip are au- 
tonomous (which we assume henceforth), i.e. for all i = 0, 1, . . . , d, VI = Vi(y) only, 
then the stochastic Taylor series for the flow is 

00 

= ^ J ai ■■ -a m (t) V ai ■■ -V am . 

m=0 a£P m 

Here P m is the set of all combinations of multi- indices a = (a±, . . . , a m ) of length m 
with oti £ {0, 1, . . . , d} and 

j ai ... am (t) = f ... f m 1 dW?* ■ ■ ■ dW^ 

Jo Jo 

are multiple Stratonovich integrals. 
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The logarithm of tp t is the exponential Lie series, Magnus expansion (Magnus 34J) 
or Chen- Strichartz formula (Chen [S], Strichartz [45 ). In other words we can express 
the flow map in the form ip t = exp ipt j where 

d d 

^ t = J2MW+ E hi J a- J Ji)(f)Wi,v j ] + --- 

i=0 j>i=0 

is the exponential Lie series for our system, and [• , •] is the Lie-Jacobi bracket on 
X(Ai). See Yamato [H], Kunita [59], Ben Arous [T] and Castell [7] for the derivation 
and convergence of the exponential Lie series expansion in the stochastic context; 
Strichartz [45] for the full explicit expansion; Sussmann [46 for a related product 
expansion and Lyons [32] for extensions to rough paths. 
Let us denote the truncated exponential Lie series by 

Ipt = 22 JaC a , (7.1) 

where Q m denotes the finite set of multi-indices a for which || J q ||l 2 is of order up to 
and including t rn , where m = 1/2, 1,3/2,.... The terms c a are linear combinations 
of finitely many (length a) products of the smooth vector fields V{, i = 0,1, ... ,d. 
The following asymptotic convergence result can be established along the lines of the 
proof for linear stochastic differential equations in Lord, Malham and Wiese [3T]; we 
provide a proof in Appendix [X] 

Theorem 7.1. Assume the vector fields Vi have 2m +1 uniformly bounded deriva- 
tives, for all i — 0,1, ... ,d. Then for t < 1, the flow expip t o yo is square-integrable, 
where ipt is the truncated Lie series (|7.1[) . Further, if y is the solution of the stochastic 
differential equation (|1.1[) . there exists a constant C(m, \\y0W2) such that 

\\y t - exp^ t o y \\ L2 < C(m, \\y \\ 2 ) t m+1 ' 2 . (7.2) 



8. Geometric Castell— Gaines methods. Consider the truncated exponential 
Lie series ipt n ,t n+1 across the interval [t n , t n+ {\. We approximate higher order multiple 
Stratonovich integrals across each time-step by their expectations conditioned on the 
increments of the Wiener processes on suitable subdivisions (Gaines and Lyons [20 ). 
An approximation to the solution of the stochastic differential equation (jl.ip across 
the interval [t n , i„+i] is given by the flow generated by the truncated and conditioned 
exponential Lie series ipt n ,t n+1 via 

Vt n+ i ~ ex P(^t„,t„+i) °Vt n ■ 

Hence the solution to the stochastic differential equation (jl.ip can be approximately 
computed by solving the ordinary differential system (see Castell and Gaines [S]; 
Misawa [Mj) 

u'(t) = ip tn ,t n+1 o u(r) (8.1) 

across the interval r S [0,1]. Then if u(Q) = yt„ we will get u(l) w Vt n+1 - We 
must choose a sufficiently accurate ordinary differential integrator to solve l|8.ip — we 
implicitly assume this henceforth. 
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The set of governing vector fields Vi, i = 0,1, ... ,d, prescribes a map from the 
driving path process w = (W 1 , . . . , W d ) to the unique solution process y £ M. to the 
stochastic differential equation (jl.ip . The map w i-^ y is called the ltd map. Recall 
that we assume the vector fields are smooth. When there is only one driving Wiener 
process (d — 1) the ltd map is continuous in the topology of uniform convergence 
(Theorem 1.1.1. in Lyons and Qian 33J). When there are two or more driving pro- 
cesses (d > 2) the Universal Limit Theorem (Theorem 6.2.2. in Lyons and Qian [33] ) 
tells us that the Ito map is continuous in the p- variation topology, in particular for 
2 < p < 3. A Wiener path with d > 2 has p- variation with p > 2, and the p- variation 
metric in this case includes information about the Levy chordal areas of the path 
(Lyons [35]). Hence we must choose suitable piecewise smooth approximations to the 
driving path process w. The following result follows from the corresponding result 
for ordinary differential equations in Hairer, Lubich and Wanner [55] (p. 112) as well 
as directly from Chapter VIII in Malliavin [36 on the Transfer Principle (see also 
Emery [15]). 

Lemma 8.1. A necessary and sufficient condition for the solution to the stochastic 
differential equation (11. ip to evolve on a smooth n-dimensional submanifold M. of 
M. N (n < N) up to a stopping time T» is that Vi(y,t) £ T y Ai for all y £ M, i = 
0,1,..., d. 

Hence the stochastic Taylor expansion for the flow ip t is a diffeomorphism on A4 . 
However a truncated version of the stochastic Taylor expansion for the flow ipt will not 
in general keep you on the manifold, i.e. if ya £ M then (p t yo need not necessarily 
lie in M.. On the other hand, the exponential Lie series ipt, or any truncation ip t of 
it, lies in X(M). By Lemma [8.11 this is a necessary and sufficient condition for the 
corresponding flow-map exp^t to be a diffeomorphism on M.. Hence if w(0) = yt n £ 
M, then yt n+1 ~ u(l) £ M. When solving the ordinary differential equation (|8.1|) . 
classical geometric integration methods, for example Lie group integrators such as 
Runge-Kutta Munthe-Kaas methods, over the interval r <E [0, 1] will numerically 
ensure yt„ +1 stays in M.. Additionally, as the following result reveals, numerical 
methods constructed using the Castell-Gaines Lie series approach can also be more 
accurate (a proof is provided in Appendix [B]) . We define the strong global error at 
time T associated with an approximate solution y-p as £ = ||yr — j/t||l 2 - 

Theorem 8.2. In the case of two independent Wiener processes and under the 
assumptions of Theorem \ 7.1\ for any initial condition yo £ M. and a sufficiently 
small fixed stepsize h = t n +\ — t n , the order 1/2 Lie series integrator is globally more 
accurate in L 2 than the order 1/2 stochastic Taylor integrator. In addition, in the 
case of one Wiener process, the order 1 and 3/2 uniformly accurate exponential Lie 
series integrators generated by iptl l ]t n+1 — JqVq + J\V\ + ^([Vi, [V\, Vq]]) and 

= JqV ° + JiVi + ^ Jni - Ji °)[ y °' ^] + t^([ Fi ' ' 

respectively, are globally more accurate in L 2 than their corresponding stochastic Tay- 
lor integrators. In other words, if £^ denotes the global error of the exponential Lie 
series integrators of order m above, and £^ is the global error of the stochastic Taylor 
integrators of the corresponding order, then £^ < £^ for m — 1/2, 1, 3/2. 
Remarks. 

1. The result for -0( 3 / 2 ) is new. That the order- 1/2 Lie series integrator (for two 
Wiener processes) and the order 1 integrator generated by t/^ 1 ' are uniformly more 
accurate confirms the asymptotically efficient properties of these schemes proved by 
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Castell and Gaines [8]. The proof follows along the lines of an analogous result for 
linear stochastic systems considered in Lord, Malham and Wiese [31] . 

2. Consider the order 1/2 exponential Lie series with no vector field commu- 
tations. Solving the ordinary differential equation (|8.ip using an (ordinary) Euler 
Munthe-Kaas method and approximating dexp" 1 w id is equivalent to the order 1/2 
stochastic Taylor Munthe-Kaas method (for the same Lie group and action). 

9. Numerical examples. 

9.1. Rigid body. We consider the dynamics of a rigid body such as a satellite 
(see Marsden and Ratiu [37] )• We will suppose that the rigid body is perturbed by two 
independent multiplicative stochastic processes W 1 and W 2 with the corresponding 
vector fields Vi(y) = £i(y) y, for i = 0, 1, 2, with £i 6 so(3). If we normalize the initial 
data yo so that \yo\ = 1 then the dynamics evolves on M. = § 2 . We naturally suppose 
Q = 50(3), and A Vo (S) = Syo so that X^(y) = £,i(y)y, and we can pull back the 
flow generated by V on M to the flow on Q generated by X^(S,t) = £ ii (^A yo (S)) S, 
i = 0,l, 2. We use the following matrix representation for the £i(y) G so(3): 

/ ~2/3M,3 2/2M,2 \ 

&(y) = V3/ai,3 -yi/ui,i , 
\-2/2M,2 yi/a M / 

where the constants cti.j for j — 1,2,3 are chosen so that the vector fields Vi and 
matrices £j do not commute for i — 0, 1,2: ao.i = 3, ao.2 = 1, ao.3 = 2, a^i = 1, 
ai ; 2 = 1/2, ai,3 = 3/2, ct2,i = 1/4, 02,2 = 1, «2,3 = 1/2. The vector fields Vi satisfy 
the conditions of Theorem 17.11 since the manifold is compact in this case. 

We will numerically solve (|1.1[) using three different order 1 methods: stochastic 
Taylor, stochastic Taylor Munthe-Kaas based on (16. 1|) and Castell-Gaines (a stan- 
dard non-geometric Runge-Kutta method is used to solve the ordinary differential 
equation (|8.1|) ). The vector field compositions ViVj needed for the stochastic Taylor 
and Castell-Gaines methods are readily computed. For the Munthe-Kaas method we 
note that we have ^00 = £,i{yo) and 

= A.{Vo,yo;oLi,aj) - !tei(yo),£j(2/o)] • 
Here o € so (3) is the zero element on the Lie algebra, and for all y, z S K 3 we define 

£i)J_/ 

and * : R 3 ^so(3) denotes the vector space isomorphism a 1— > <r where 

(0 -CT3 (72 \ 
cr 3 -(Ji . 
-(72 (Ti / 

Note that yz = y A z (see Marsden and Ratiu [37] ). Note also since <r 6 so (3), 
exper G 50(3) can be conveniently and cheaply computed using Rodrigues' formula 
(see Marsden and Ratiu [37] or Iserles et al. |25j). 

In Figure I9TTI we show the distance from the manifold § 2 of each the three approx- 
imations; we start with initial data yo — (\/2, a/2, 0) t . The stochastic Taylor Munthe- 
Kaas method can be seen to preserve the solution in the unit sphere to within machine 
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Fig. 9.1. Rigid body: We show the log-distance of the approximate solution to the unit sphere 
as a function of time for each of the methods. Below we show the approximate solutions as a function 
of time for the stochastic Taylor (blue) and Munthe-Kaas methods (magenta). The trajectory starts 
at the top right and eventually drifts over the left horizon. 
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error. We also see that the stochastic Taylor method clearly drifts off the sphere as 
the integration time progresses, as does the non-geometric Castell-Gaines method — 
which does however remain markedly closer to the manifold than the stochastic Taylor 
scheme. 
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Fig. 9.2. Autonomous underwater vehicle: We show the log-distance of the approximate solution 
to the two Casimirs Ci = it ■ p ( dotted line) and C2 = \p\ 2 (solid line) as a function of time for each 
of the methods. Below, we also show the global error as a function of stepsize. 
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9.2. Autonomous underwater vehicle. The dynamics of an ellipsoidal au- 
tonomous underwater vehicle is prescribed by the state y = (tt,p) G se(3)* where 
7r G so (3)* is its angular momentum and p G (R 3 )* its linear momentum (see Holmes, 
Jenkins and Leonard [24], Egeland, Dalsmo and S0rdalen [12] and Marsden and 
Ratiu [37)). We suppose that the vehicle is perturbed by two independent multi- 
plicative stochastic processes. The governing vector fields are for i = 0, 1, 2: 

Vi(y) = ad*, oy. 

Here = G se(3) where a^(y) = /r 1 ^ and = iWT x p are the 

angular and linear velocity, and Ii — diag(ai i i, a^2, 0^,3) and Mi = diag(/3i i, /?i,2, ^,3) 
are the constant moment of inertia and mass matrices, respectively. Explicitly for 
£ G se(3) we have 

ad^ o y = (it A lu + p A u, p A u>) . 

The system of vector fields Vi, i — 0,1,2 represents the Lie-Poisson dynamics on 
M. = se(3)* (Marsden and Ratiu [37]). There are two independent Casimir functions 
Cfe : se(3)*— *E, k = 1,2, namely C\ = it ■ p and C 2 = \p\ 2 ; these are conserved by the 
flow on se(3)*. Note that the Hamiltonian, i.e. total kinetic energy ^(tt ■ uj + p ■ u), 
is also exactly conserved (and helpful for establishing the sufficiency conditions in 
Theorem 17. ip . but that is not our focus here. 

If G = SE(3) ~ SO(3) x M 3 , then the coadjoint action of SE(3) on se(3)*, 
Ad* : SE(3) x se(3)*^se(3)* is defined for all S = (s,p) G SE(3), where s G SO(3) 
and p G M 3 , and y G se(3)* by: A y o S — AdJ-i o y = (s7r + p A (sp),sp). The 
corresponding infinitesimal action A: se(3) x se(3)*^se(3)* for all £ G se(3) and y G 
se(3)* is given by (see Marsden and Ratiu [37], p. 477) 

h y = ~ad*- o y ■ 

Since ad*-(y) = — Ag(j/) = X-^(y) the governing set of vector fields on se(3)* are 

Vi(y) = A_ Si oy. 

We can now pull back this flow on se(3)* to a flow on SE(3) via Aj, . The correspond- 
ing flow on SE(3) is generated by the governing set of vector fields for i = 0, 1, 2: 

oS = -(u)i(y) A s,uji(y) Ap + u,(j/)) , 

with y = A yo (5'). 

To aid implementation note that SE(3) = {(s,p) G SE(3): s G SO(3), p G M 3 } 
embeds into S'L(4;]R) via the map 

S = (S,p) ' ^ ^ , 

where O is the three- vector of zeros. Also se(3) is isomorphic to a Lie subalgebra of 
s[(4;R) with elements of the form 
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Hence the governing vector fields on SE(3) are of the form = — 5, where 

The governing vector fields on se(3) are v^(a) = — dexp CT o £j (A^ (exp a)) . Again the 
vector field compositions VtVj needed for the stochastic Taylor and Castell-Gaines 
methods can be computed straightforwardly. Direct calculation also reveals that in 
block matrix form 

v^oo = ( A ^> *<* ai > <*i)+ T A(po,Po; ft, «i) A(7r ,po; a u ^ 

Here A(y, z; a, j3) is defined as for the rigid body example. Note that the exponential 
map exp se ( 3 ) : se(3)^SE(3) is defined for all a = (9, £) e se(3) by 

where exp sc ,( 3 ) is the exponential map from so (3) to 5*0(3) which can be computed 
using Rodrigues' formula and (see Bullo and Murray 0], p. 5) 

/(*) = h*s + (i - cos \\e\\)§/\\er + (i - (sin \\e\\)/\\e\\)ey\\e\\ 2 . 

In Figure Wl2\ we show the distance from the manifold se(3)* of each the three ap- 
proximations; in particular how far the individual trajectories stray from the Casimirs 
C\ = 7r ■ p and C 2 = |j>| 2 - We start with the initial data y = (y/2, y/2, 0, 0, y/2, V%) T ■ 
As before the stochastic Taylor Munthe-Kaas method can be seen to preserve the 
Casimirs to within machine error. We also see that the stochastic Taylor method 
clearly drifts off the manifold as the integration time progresses and at a particular 
time depending on the Wiener path shoots off very rapidly away from the manifold. 
Note also that for large stepsizes the stochastic Taylor method is unstable. However 
the non-geometric Castell-Gaines and stochastic Munthe-Kaas methods still give reli- 
able results in that regime. Lastly, although the the stochastic Munthe-Kaas method 
adheres to the manifold to within machine error, the error of the non-geometric 
Castell-Gaines method is actually smaller. 

10. Conclusions. We have established and implemented stochastic Lie group 
integrators based on stochastic Munthe-Kaas methods and also derived geometric 
Castell-Gaines methods. We have also revealed several aspects of these integrators 
that require further investigation. 

1. We could construct a stochastic nonlinear Magnus method by approximating 
the solution to the stochastic differential equation (|5.ip on the Lie algebra using Picard 
iterations (see Casas and Iserles [6]). 

2. We would like to develop a practical procedure for implementing ordinary 
Munthe-Kaas methods for higher order Castell-Gaines integrators. We need to de- 
termine the element £ : M. — >g so that in (|8.ip we have i/j — A £ . 

3. We need to determine the properties of the local and global errors for the 
stochastic Munthe-Kaas methods. Also a thorough investigation of the stability prop- 
erties of the stochastic Munthe-Kaas and Castell-Gaines methods is required. For 
the autonomous underwater vehicle simulations they were both superior to the direct 
stochastic Taylor method, especially for larger stepsizes. We also need to compare the 
relative efficiency of the methods concerned, in particular to compare an optimally 
efficient geometric Castell-Gaines method with the stochastic Munthe-Kaas method. 
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4. Although we have chiefly confined ourselves to driving paths that are Wiener 
processes, we can extend Munthe-Kaas and Castell-Gaines methods to rougher driv- 
ing paths (Lyons and Qian [33], Friz [TS], Friz and Victoir [H]). Further, what hap- 
pens when we consider processes involving jumps? For example Srivastava, Miller and 
Grenander [44] consider jump diffusion processes on matrix Lie groups for Bayesian 
inference. Or what if we consider fractional Brownian driving paths; Baudoin and 
Coutin [3J investigate fractional Brownian motions on Lie groups? 

5. Schiff and Shnider 43] have used Lie group methods to derive Mobius schemes 
for numerically integrating deterministic Riccati systems beyond finite time removable 
singularities and numerical instabilities. They integrate a linear system of equations 
on the general linear group GL(n) which corresponds to a Riccati flow on the Grass- 
mannian manifold Gr(fc, n) via the Mobius action map. Lord, Malham and Wiese [31j 
implemented stochastic Mobius schemes and show that they can be more accurate and 
cost effective than directly solving stochastic Riccati systems using stochastic Taylor 
methods. We would like to investigate further their effectiveness for stochastic Ric- 
cati equations arising in Kalman filtering (Kloeden and Platen [27] ) and to backward 
stochastic Riccati equations arising in optimal stochastic linear-quadratic control (see 
for example Kohlmann and Tang [28] and Estrade and Pontier [16]). 

6. Other areas of potential application of the methods we have presented in this 
paper are for example: term-structure interest rate models evolving on finite dimen- 
sional invariant manifolds (see Filipovic and Teichmann [17]); stochastic dynamics 
triggered by DNA damage (Chickarmane, Ray, Sauro and Nadim [10]) and stochastic 
symplectic integrators for which the gradient of the solution evolves on the symplectic 
Lie group (see Milstein and Tretyakov [55]). 
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Appendix A. Proof of Theorem (77TJ We follow the proof for linear stochastic 
differential equations in Lord, Malham and Wiese [31] (where further technical details 
on estimates for multiple Stratonovich integrals can be found). Suppose ip t = ip t {m) 
is the truncated Lie series (]7.1|) . First we show that exp ipt ° yo G L 2 . We see that 
for any number k, [ipt) ° Ho is a sum of |Q m | fe terms, each of which is a fc-multiple 
product of terms J a c a o y . It follows that 



Note that the maximum of the norm of the compositions of vector fields c a oy is taken 
over a finite set. Repeated application of the product rule reveals that for i = 1, . . . , k, 
each term ' J ai J a2 ■ ■ ■ J ak ' in (jA.l[) is the sum of at most 2 2mk ~ 1 Stratonovich integrals 
J/3, where for t < 1, \\Jp\\ L i < 2 4mfe ~ 1 t k / 2 . Since the right hand side of equation (jA~T|) 
consists of l^ml* 2 2mk ~ 1 Stratonovich integrals Jp, we conclude that, 




(A.l) 



--l,...,k 



(A) k oy <( max \\c a o y || • |Q m | . 2 6m • t 1 

L 2 \o6(J m 



) 



k 
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Hence exp ip t o y is square-integrable. 

Second we prove (|7.2p . Let y t denote the stochastic Taylor series solution, trun- 
cated to included terms of order up to and including t m . We have 

\\y t - exp^ o y \\ L2 < \\y t - y t \\ L2 + \\y t - exp ^ o y \\ L2 . 

We know y t € L 2 — see Lemma III. 2.1 in Gihman and Skorohod [21]. Note that 
the assumptions there are fulfilled, since the uniform boundedness of the derivatives 
implies uniform Lipschitz continuity of the vector fields by the mean value theorem, 
and uniform Lipschitz continuity in turn implies a linear growth condition for the 
vector fields since they are autonomous. Note that yt is a strong approximation to 
yt up to and including terms of order t m , with the remainder consisting of 0(t m+1 ^ 2 ) 
terms (see Proposition 5.9.1 in Kloeden and Platen [27] ). It follows from the definition 
of the exponential Lie series as the logarithm of the stochastic Taylor series, that the 
terms of order up to and including t m in exp^t o y correspond with yt; the error 
consists of 0{t m+1 / 2 ) terms. 

Appendix B. Proof of Theorem 18.21 Our proof follows along the lines of that 
for uniformly accurate Magnus integrators for linear constant coefficient systems (see 
Lord, Malham & Wiese [3T] and Malham and Wiese [35]). Let tpt n ,t„ +1 ar >d ( fit n ,t rl+1 
denote the exact and approximate flow-maps constructed on the interval [t„, t n+ {\ of 
length h. We define the local flow remainder as 



+i ' 



and so the local remainder is Rt n ,t„ +1 yt n - Let R ls and R st denote the local flow 
remainders corresponding to the exponential Lie series and stochastic Taylor approx- 
imations, respectively. 

B.l. Order 1/2 integrator: two Wiener processes. For the global order 1/2 
integrators we have to leading order i? ls = \{ J\i — J2i)\V\, V2] and R st — J\iV\V2 + 
■/21V2V1. Note that we have included the terms J\\V 2 and J22V2 in the integrators. 
A direct calculation reveals that 

E((i? st o y ) T R st o y ) = E((i? ls o y ) T R ls o y ) + h 2m U T BU + 0(h 2m+ ^) . (B.l) 

Here m = 1/2 (for the order 1/2 integrators), U = {ViV 2 o y , V 2 Vi o y Q ) T e R 2n , and 
B e Jj2nx2n cons j s t s of H x ii diagonal blocks of the form hjlnxn where 

b=i( 1 1 
* \1 1 

and I n xn is the nxn identity matrix. Since b is positive semi-definite, the matrix B = 
b®I n xn is positive semi-dcfinitc. Hence the order 1/2 exponential Lie series integrator 
is locally more accurate than the corresponding stochastic Taylor integrator. 

B.2. Order 1 integrator: one Wiener process. For the global order 1 in- 
tegrators we have to leading order R ls = |(Joi — ^io)^, ^1] an d R st = JqxVqVi + 
JwV\V + JinV^ + \h 2 (V V 2 + V 2 V ). The terms of order h 2 shown are significant 
when we consider the global error in Section IB. 41 below. The estimate (IB.1[) also 
applies in this case with m = 1 and U — (V0V1 o y , V1V0 o y , o y ) T g R 3n ; and 
B e Jj3nx3n cons i s ts of n x n diagonal blocks of the form bijl nxn where 
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Since b is positive semi-definite, the matrix B = b® I n xn is positive semi-definite. 
Hence the order 1 exponential Lie series integrator is locally more accurate than the 
corresponding stochastic Taylor integrator. 

B.3. Order 3/2 integrator: one Wiener process. The local flow remainders 
are R ls = \ ( J 110 -2 J wl + J 011 - \h 2 ) [V u [V u V }] and R st = J 011 V Q V 2 + J 10 iVi Wi + 
•hwVfVo + J nll Vf - \h 2 (V V^ + V?V + \V?). The terms of order h 2 shown are 
significant when we consider the global error — but for a different reason this time — see 
Section [B.4I below. Again, the estimate (|B.1[) applies in this case with to = 3/2 and 
U = (V V? o y , ViV Vi o j/o, V?V o z/o, Vj 4 o y ) T e K 4n ; and B e R 4 « x4 « consists of 
n x n diagonal blocks of the form bijI n - Kn where 

^11 8 5 12^ 
8 8 8 12 
5 8 11 12 
\12 12 12 24/ 

Again, B is positive semi-definite and the order 3/2 exponential Lie series integrator 
is locally more accurate than the corresponding stochastic Taylor integrator. 

B.4. Global error. Recall that we define the strong global error at time T 
associated with an approximate solution yx as £ = \\yr — Vt\\l 2 - The exact and 
approximate solutions can be constructed by successively applying the exact and 
approximate flow maps ift n ,t n+1 and 0t„,t n+1 on the successive intervals [i„,i„+i] to 
the initial data ?/o- A straightforward calculation shows for a small fixed stepsize h, 

£ 2 =E{Koy Q ) T Koy Q: (B.2) 

up to higher order terms, where TZ = ^2^=0 Pt n +i,tN ° Rt n ,t n +i ° fto.tn is the standard 
accumulated local error contribution to the global error. The important conclusion is 
that when we construct the global error (|B.2|) , the terms of leading order in the local 
flow remainders R ls or R st with zero expectation lose only a half order of convergence 
in this accumulation effect. Hence in the local flow remainders shown above, for the 
terms of zero expectation, the local superior accuracy for the Lie series integrators 
transfers to the corresponding global errors (see Lord, Malham and Wiese [31] for 
more details). Terms of non-zero expectation however behave like deterministic er- 
ror terms losing a whole order (in the local to global convergence); they contribute 
to the global error through their expectations. Hence we include such terms of or- 
der h 2 in the order 3/2 integrators above and they appear as the terms subtracted 
from the remainders shown. For the order 1 integrators we do not need to include 
the order h 2 terms in the integrator to obtain the correct mean-square convergence. 
However to guarantee that the global error for the exponential Lie series integrator 
is always smaller than that for the stochastic Taylor scheme, we include this term in 
the integrator. 
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